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13  ABSTRACT 


Dual  variational  principles  for  steady  state  wave  propagation  in  three  dimen¬ 
sional  thermoviscoelastic  media  are  presented.  The  first  one,  for  the  equations  of 
motion,  involves  only  the  complex  displacement  function.  The  second  principle  is 
for  the  energy  equation.  The  specialized  versions  of  these  principles  in  two- 
dimensional  polar  coordinates  and  then  in  one  dimension  are  obtained.  A  one-dimen¬ 
sional  example,  that  of  wave  propagation  in  a  thermoviscoelastic  rod  insulated  on 
its  lateral  surface  and  driven  by  a  sinusoidal  stress  at  one  end,  is  solved  using 
the  Rayleigh-Ritz  method.  The  displacement  and  temperature  functions  are  expressed 
as  series  of  polynomials.  Successive  approximations  for  the  solution  are  compared 
with  a  solution  obtained  by  a  method  of  finite  differences,  and  an  estimate  of  the 
degree  of  accuracy  as  a  function  of  the  number  of  terms  taken  in  the  series  is  ob¬ 
tained.  It^ is  found  that  as  long  as  the  spatial  distribution  of  stress  and  tempera¬ 
ture  are  sufficiently  smooth,  r  ...  convergence  to  the  correct  solution  is  obtained. 
If  the  stress  is  a  rapidly  oscii.  ting  function  of  the  distance  along  the  rod,  poly¬ 
nomials  are  no  longer  efficient  and  other  test  functions  must  be  chosen. 
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INTRODUCTION 


i 


Several  authors  have  developed  and  used  variational  principles  to 

obtain  solutions  to  problems  in  quasi-static  and  dynamic  viscoelasticity, 

* 

with  and  without  thermomechanical  coupling.  Gurtin  [1]  and  Leitman  [2] 
have  developed  variational  principles  for  viscoelastic  media  without  ther¬ 
momechanical  coupling.  They  have  used  the  convolution  form  of  the  consti¬ 
tutive  equations  and  have  developed  variational  principles  for  several 
types  of  boundary  value  problems.  Their  work  appears  to  be  primarily  of 
mathematical  interest.  Valanis  [3]  has  developed  a  principle  applicable 
to  viscoelastic  materials  with  constant  Poisson's  ratio,  without  thermo¬ 
mechanical  coupling. 

Schapery  [6,7]  has  studied  wave  propagation  in  viscoelastic  media 
with  thermomechanical  coupling.  In  [7]  he  has  used  a  complex  modulus  form 
of  the  constitutive  equations  and  has  developed  a  variational  principle 
analogous  to  Reissner's  complementary  principle  using  complex  kinetic  and 
potential  "energy"  functions.  His  principle,  however,  involves  both  stress 
and  displacement  functions  which  must  already  satisfy  the  equations  of  mo¬ 
tion.  He  has  considered  examples  with  bodies  that  are  either  massless  or 
with  concentrated  mass,  and  in  his  last  example  of  a  'solid  cylinder  with 
distributed  mass'  he  only  gives  a  first  approximation  to  the  solution, 
using  only  one  term  of  a  series  expansion.  While  his  method  appears  promis¬ 
ing,  the  question  of  convergence  to  the  exact  solution,  or,  in  other  words, 
how  many  terms  in  the  series  are  necessary  to  get  a  sufficiently  close 
approximation  to  the  exact  solution,  remains  open. 

* 

Numbers  in  square  brackets  designate  references  at  the  end  of  this 
report. 
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Tins  report  is  concerned  with  the  application  of  variational  prin¬ 
ciples  to  problems  of  steady  state  wave  propagation  in  viscoelastic  media 
with  thormomoclianical  coupling.  A  complex  modulus  description  of  the 
constitutive  equations  is  used.  The  material  is  assumed  to  be  thermorheo- 
logically  simple  f 5 ]  and  the  energy  equation,  as  suggested  by  Schapery 
1 7  I .  uses  the  cycle  averaged  temperature  distribution  with  the  cycle 
averaged  value  of  the  Rayleigh  dissipation  function  acting  as  the  heat 
source.  The  displacement  variational  principle  suggested  here  involves 
only  the  complex  displacement  function  and  an  admissible  set  of  dis¬ 
placement  functions  need  only  satisfy  any  prescribed  displacement  boundary 
condition:;  LL.ic  might  exist.  This  principle  can  be  considered  to  be  an 
extension  of  that  developed  by  Kohn,  Krumhnnsl  and  Lee  f4]  for  elastii 
media.  Tt  uses  complex  instead  of  real  "energy"  functions.  The  tem¬ 
perature  variational  principle  is  the  one  suggested  by  Biot  F 8 ]  and 
Schapery  f  7 1 . 

An  alternative  form  of  these  principles  is  suggested.  This  proves 
more  useful  for  certain  applications.  These  principles  are  set  up  for 
general  three-dimensional  problems  and  are  later  specialized  to  the 
cases  of  two  and  one  dimension. 

As  an  example,  the  problem  of  steady  state  longitudinal  waves  in  a 
viscoelastic  rod  with  thermal  coupling  subjected  to  a  sinusoidal  stress 
applied  at  one  end,  is  solved  using  a  variational  approach.  Huang  and 
Leo  19]  solved  this  problem  including  time  as  an  independent  variable. 

This  resulted  in  partial  differential  equations  which  were  solved 
numerically  using  a  method  of  finite  differences.  This  is  useful  if 
the  time  histories  of  the  stress  and  temperature  have  to  be  determined. 
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For  most  engineering  design  applications,  however,  the  steady  state 
values  of  stress  and  temperature  are  of  primary  interest,  since  due  to 
dissipation  of  mechanical  energy  the  temperature  increases  until  a  steady 
state  is  reached,  if  in  fact  the  situation  is  stable.  Such  a  steady 
i  state  yields  the  most  severe  temperature  conditions  which  are  the  major 

r 

i 

concern  in  design.  In  such  cases  it  is  simpler  and  far  more  efficient 
\  to  obtain  the  steady  state  values  directly  instead  of  following  the  com¬ 

plete  time  history  of  the  process  till  a  steady  state  is  reached.  In 
this  example,  the  steady  state  values  of  stress  and  temperature  have 
been  directly  obtained  by  using  a  Rayleigh-Ritz  procedure  on  the  alter¬ 
native  form  of  the  variational  principles.  Functions  for  displacement 
,  (complex)  and  temperature  are  assumed  as  polynomial  series  (for  con¬ 

venience)  with  "n"  and  "m"  terms  respectively,  with  unknown  coefficients. 
Simultaneous  extremization  of  two  functionals  is  carried  out  by  solving 
the  resultant  nonlinear  algebraic  equations  in  a  computer.  The  number 
of  terms  "n"  and  "m"  can  be  set  in  the  program.  Calculations  for  a 
Lockheed  solid  propellant  [9]  are  carried  out  for  various  values  for 
"n"  and  "m"  and  the  question  of  rapidity  of  convergence  to  the  solution 
given  in  [9]  is  discussed. 


GOVERNING  DIFFERENTIAL  EQUATIONS 


1.  General  Equations  in  Three  Dimensions 

Let  us  consider  the  governing  differential  equations  for  stresses, 
displacements  and  temperature  in  steady  state  oscillations  of  linear 
isotropic  viscoelastic  media.  The  thermomechanical  coupling  is  caused 
by  the  cycle  averaged  value  of  the  mechanical  dissipation  function  acting 
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as  the  heat  source  it;  J  -s  energy  equation  and  by  the  fact  that  the  com¬ 
plex  viscoelastic  modulii  are  temperature  dependent.  As  pointed  out  by 
Schapery  [6,7]  the  coupling  terms  due  to  the  dilatation  and  potential 
energy  drop  out  of  the  energy  equation  if  it  is  integrated  over  a  cycle. 

x,'!'-  assume  steady  sta  a  conditions  where  the  mechanical  variables 
are  <•.•  -.’ionic  fund  Ions  of  time,  and  the  temperature,  after  a  sufficiently 
long  ? me,  is  i.vtepeodont  of  time.  Strictly  speaking,  the  temperature  is 
never  truly  tim^  ].».  dependent  but  has  small  cyclic  variations  about  a  mean 
value  as  a  result  of  the  cyclic  variations  of  the  potential  energy, 
dilatation  and  dissipation  (see  [9]).  These  small  fluctuations,  however, 
will  be  neglected  and  henceforth  the  temperature  will  mean  its  cycle 
averaged  steady  state  value. 

Let  the  stress  and  strain  tensors  and  the  displacement  vector  be 
defined  as  the  real  parts  of 
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where  i  =  u>  is  the  frequency  (real)  and  t  is  time.  The 

complex  quantities  and  u^  are  used  most  of  the  time  in 


further  calculations  and  will  be  referred  to  as  simply  stress,  strain 
and  displacement  respectively.  The  familiar  cartesian  tensor  notation 
is  used  here.  The  suffixes  i  and  j  range  from  1  to  3  and  a  repeated 
index  implies  summation  over  that  index. 

The  equations  of  motion  are 


ZiiSlU 


.  :rni- 


-rrirlUr’i  n  ■»-  -  ?ww.>  ^ 


3* 


i 


i 


3 


Jm 


:wstwrw«:^»3®«ik\  *»  vwi. 


aij,5  +  p  "  ui  *  0 


where  P  is  the  mass  density  (assumed  constant). 

The  constitutive  equations,  using  the  familiar  complex  modulus 


formulation,  are 


“ij'*-  "k,k6ij  +  i"  (“i,j  +  uj,i>  <3) 

where  \  and  p.  are  complex  Lame  parameters  which  are  functions  of 
temperature,  and  is  the  Kronicker  delta. 

Equation  (3)  is  entirely  analogous  to  the  constitutive  equations 
for  elastic  solids  and  is  obtained  directly  by  the  use  of  the  well 
known  correspondence  principle  [12]. 

We  note  that  Eqs.  (3)  should  contain  a  term  due  to  thermal  expan¬ 
sion.  However,  as  pointed  out  by  Schapery  in  [7],  in  our  problem  each 
of  the  mechanical  variables  can  be  separated  into  two  parts:  (a)  that 
due  to  applied  cyclic  loading  in  the  absence  of  thermal  expansion,  and 
(b)  that  due  to  thermal  expansion  with  homogeneous  mechanical  boundary 
conditions  and  a  temperature  distribution  obtained  from  the  components 
of  (a) .  Here  we  are  only  concerned  with  the  former  components  and  the 
latter  can  be  obtained  from  standard  thermoelastic  analysis. 

The  Lame  parameters  are  related  to  the  more  commonly  used  complex 
shear  and  bulk  compliances  J  and  B  by  the  relations 
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Typically^  in  polymers,  J  is  a  very  strong  nonlinear  function 

•ft 

of  temperature  while  B  is  a  relatively  weak  function  of  temperature. 

* 

In  thermo rheological ly  simple  materials  it  is  assumed  that  J  is 
a  function  of  only  the  reduced  frequency  <i>*  which  is  related  to  the 
actual  frequency  u>  through  the  temperature  dependent  shift  factor 
aT  (see  [7]),  i.e. 


u)'  =  u)  a^T) 


(5) 


where  a^,  represents  the  effect  of  temperature  on  viscosity. 

Combining  Eqs.  (2)  and  (3)  one  can  write  the  equations  of  motion  in 
terms  of  displacement  alone 


<’*  “k.k’.i  + 


The  steady  state  energy  equation  for  the  cycle  averaged  temperature 
distribution  is  given  by 


KT,ii  =  -2D  (7) 

where  T  is  the  temperature,  K  is  the  thermal  conductivity  (assumed 
constant)  and  D  is  the  cycle  averaged  value  of  the  Rayleigh  dissipa¬ 
tion  function  given  by 

dt» 

where  Re  denotes  the  real  part  of  the  complex  argument.  2D  is  the 
cycle  averaged  value  of  the  mechanical  dissipation. 


Using  Eqs.  (1),  (3)  and  (7)  and  after  carrying  out  the  necessary 
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integration,  we  have 


B-SW\ki2+*2|eijl2j 


where 


X  =  X^  +  iX^ 


ii  «=  |X1  4  iM-2 


X^  and  ^  are  the  storage  Lame  parameters  and  X2  ,  H2  the  loss 


parameters . 
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and  "— "  denotes  the  complex  conjugate. 

Note  that  as  defined  D  is  a  real  function  of  the  strain  tensor 
and  the  loss  parameters  X2  and  |i2  . 

Since  X2  and  p,2  are  nonlinear  functions  of  temperature,  (7)  is 
a  nonlinear  partial  differential  equation  for  the  temperature.  Also, 
the  steady  state  temperature  is  purely  a  spatial  function  independent 
of  time. 

Using  the  familiar  kinematic  relations 


€1J  “  2<ui,j  +  “j,i>  <9) 

we  can  write  Eq.  (7)  in  terms  of  displacement  and  temperature 

K  T,u  +!  lX2 K.k‘t 2  +  »2(ul.J  “l.J  +  ul,j  ■  0  <10> 

Equations  (6)  and  (10)  are  a  complete  set  of  four  nonlinear  partial 
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differential  equations  for  the  four  unknowns  u^  (1*1,3)  and  T  . 
Since  these  are  written  in  terms  of  displacements,  the  compatibility 
conditions  are  automatically  satisfied. 

For  the  boundary  conditions  we  assume  that  the  displacement  vector 
is  prescribed  on  a  part  of  the  Surface  Au  ,  while  the  traction  vector 
is  prescribed  on  the  remainder  Aa  .  Also,  the  temperature  is  pre¬ 
scribed  on  the  portion  of  the  surface  A^,  and  the  heat  flux  (per  unit 
area)  is  prescribed  on  the  remaining  surface  A„  . 


2.  Equations  in  Two-Dimensional  Polar  Coordinates 

It  is  instructive  to  look  at  the  form  of  the  general  equations  for 
the  case  of  plane  strain  in  polar  coordinates.  These  equations  and 
the  associated  variational  principles  to  be  discussed  later  are  useful 
for  problems  involving  long  circular  cylinders  in  plane  strain. 

The  equations  of  motion  now  take  the  well  known  form 
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while  the  constitutive  euqations  become 

cr  =  **<«r  +  V  +  er 

cQ  =  \*(er  +  e0)  +  i*  €q  (12) 
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where  a^  ,  ,  a^g  are  the  stress  components,  ,  «g  t  e^g 


the 


•*  S,  e^«fiV»VTfe  m*  V>«?v>S2rv  •. 

BSyWSSWTfSE!!^^ 


*9- 


and  u  ,  u0  the  displacement  components  in  polar 


strain  components  ana  uf  »  «g 
(r,S)  coordinates. 

The  strains  and  displacements  are  related  by  the  familiar  kinematic 
equations 


du 


er  =  dr 
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(13) 
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The  energy  equation  (Eq.  (7))  takes  the  form 


K  v  T  **  -2D 


(14) 


where 


2  d2  +ld.+i^Si 

•-  +  T  dr  +  r2  ae2 


dr 


D  =  \  [(X2  +  2^2)(er  \  +  ®q  V  +  ^2(er  *9  +  *9  V 


+  ^2  6r0  ers] 


This  can  easily  be  shown  by  using  the  fact  that  |ekk  |  and 


|6i  j  2  are  invariant  under  the  transformation  of  coordinates  and 


KJ2->  (er  +  eeH7r+I9> 
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VARIATIONAL  PRINCIPLES 


1.  Variational  Principles  in  Three  Dimensions 

The  field  equations  (6)  and  (10)  for  the  coupled  thermomechanical 
problem  together  with  the  boundary  conditions  are  equivalent  to  two 
variational  principles. 

The  variational  principle  for  the  equations  of  motion  and  the  con¬ 
stitutive  equations  (2)  and  (3)  can  be  stated  as:  of  all  displacement 
functions  u^  satisfying  prescribed  displacements  u^  on  Au  ,  the 
displacement  function  satisfying  the  equations  of  motion  (2),  the 
constitutive  equations  (3)  and  the  traction  boundary  condition  on 
is  determined  by 


a 


<uv  -  K  )  iV 
u  -V  A. 


3lj  ”j 


=  0 


(15) 


where  Uv  and  are  analogous  to  the  elastic  strain  energy  density 
and  kinetic  energy  density  and  are  given  by 

lvsl"2,1"1 


© 

n^  =  prescribed  traction  on  Aa 

n,  »  direction  cosines  of  the  outward  unit  normal  to  the  surface  A_ 
J  o 

6  means  that  the  variations  must  be  taken  with  respect  to  the  dis- 
u 


placement  function  only.  • 
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If  the  kinematic  equations  (10)  are  also  satisfied  (i.e.,  we  define 
strain  functions  to  satisfy  Eq.  (10)),  we  can  write 

•k 

\  Ok 

U  =  — - Ye  )*  +  a  e  e 
v  2  *  kk'  ^  ij  ij 


This  variational  principle  is  analogous  to  Hamilton's  principle  in 
dynamic  elasticity. 

Comparing  Eq.  (15)  with  the  variational  principle  given  by  Schapery 
in  [7]  we  see  that  here  we  can  choose  trial  functions  for  the  displace¬ 
ment  which  need  only  satisfy  the  displacement  boundary  conditions  of 
the  problem  whereas  in  [7]  Schapery  must  choose  displacement  and  stress 
functions  which  must  already  satisfy  the  equations  of  motion.  The 
latter  principle  thus  appears  more  restrictive  and  would  be  more  dif¬ 
ficult  to  apply  in  complicated  problems. 

The  variational  principle  can  be  proved  by  carrying  out  variations 
with  respect  to  the  function  u^  to  yield 


/  {(>*  uk,k\i  +  (y*  ui,j>,j  +  <“*  uj,i>,j  +  p  0,2  uj  6  ui  dV 

V 

/  /X*  \,k  6ij  +  “^i, j  +  Uj,i)  •  *ij t  nj  6  Ui  ^  "  0  <16> 


In  view  of  the  arbitrariness  of  5u^  ,  this  expression  equals  zero 
only  if  the  equations  of  motion  (2),  the  constitutive  equations  (3), 
and  the  traction  boundary  conditions  on  Aa  are  satisfied. 

If  we  restrict  the  admissible  class  of  displacement  functions  such 
that  the  boundary  conditions  for  both  the  displacements  on  Au  and 
traction  on  Aa  (through  Eq.  (3))  are  satisfied,  the  surface  integral 


drops  out  of  Eq.  (16)  and  we  are  left  with  a  simplified  form  of  the 
principle 

(  1<x*  uk,k\i  +  (11*  ui,j>,j  +  (|1*  +  p  “2  “J  6  ui  dV  ■ 0 

Jv 


Equation  (16)  (or  (17)  which  is  a  special  case  of  (16))  can  be  con¬ 
sidered  to  be  an  alternative  form  of  the  variational  principle  (15). 
Equation  (17)  can  be  considered  to  be  a  Galerkin  formulation  of  the 
problem. 

It  is  useful  to  compare  the  relative  advantages  of  the  two  forms. 
Equation  (15)  uses  energy  invariants  and  therefore  appears  more  con¬ 
venient  in  complicated  coordinate  systems.  However,  when  carrying  out 
a  Rayleigh-Ritz  method  of  solution,  use  of  Eq.  (16)  can  save  a  large 
amount  of  calculations  since  the  variations  have  already  been  carried 


•jf  ^ 

It  must  be  remembered  that  X  and  n  are  temperature  dependent 
and  in  order  to  get  the  temperature  field  we  require  another  variational 
principle  from  the  energy  equation.  This  can  be  stated  as  follows: 
of  all  temperature  distributions  which  satisfy  prescribed  T  on  , 
the  temperature  distribution  which  also  satisfies  the  energy  equation 
(7)  and  the  heat  flow  boundary  condition  on  is  determined  by 


6  {  f  tST  "  V  dV  +  f  H  T  dA;  -  0  (18) 

1  Jv  \ 

where  S^,  is  proportional  to  the  entropy  production  density  resulting 
from  temperature  gradients  (see  [7]) 


ST  -  ^  KT>i  T#i 


and  SM  is  the  integral  of  the  mechanical  dissipation 


H  =  prescribed  heat  flux  per  unit  area  out  of  the  body. 


6  means  the  variations  must  be  taken  with  respect  to  the  tempera- 
T 

ture  only; 


D  is  as  given  by  Eqs.  (7)  and  (10), 


This  principle  can  be  proved  by  taking  variations  with  respect  to 


T  to  yield 


-  I  (2D  +  KT,ii>  6T  dV  +  j  (K  n±  +  H)  6t  dA  -  0 

•V  *4,. 


In  view  of  the  arbitrariness  of  6t  ,  this  expression  is  zero 
only  if  the  energy  equation  (10)  and  the  heat  flow  boundary  condition 
on  is  satisfied. 

If  we  restrict  the  admissible  class  of  temperature  functions  such 
that  the  boundary  conditions  for  both  the  temperature  on  A^,  and  heat 
flow  on  Ajj  are  satisfied,  the  surface  integral  drops  out  of  Eq.  (19) 


and  we  are  left  with 


X(2D  +  KT,ii)  6T  dV  -  0 
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Equation  (19)  is  an  alternative  form  of  the  variational  principle  (18). 
Equation  (18)  uses  thermodynamic  invariants  and  the  comments  made  ear¬ 
lier  about  the  two  forms  of  the  variational  principle  for  the  equa¬ 
tions  of  motion  apply  here  too. 

The  variational  principles  for  displacement  and  temperature  (Eqs. 
(15)  and  (18))  are  entirely  equivalent  to  the  field  equations  (2),  (3), 
and  (10),  with  their  associated  boundary  conditions.  The  displacement 
and  temperature  functions  can  be  obtained  by  simultaneously  making  the 
appropriate  integrals  stationary  with  respect  to  displacement  and  tem¬ 
perature  respectively.  The  first  equation  (15)  could  be  regarded  as 
getting  a  stationary  "cost"  function,  and  the  second  (18)  as  a  con¬ 
straint,  or  vice  versa. 

A  Rayleigh-Ritz  procedure  can  be  used  to  obtain  the  displacement 
(and  hence  the  stress  components)  and  temperature  distribution.  This 
is  done  in  a  one- dimensional  example  presented  later  in  the  report. 


2.  Variational  Principles  in  Two-Dimensional  Polar  Coordinates  for 
Plane  Strain 

The  variational  principle  for  displacement  (Eq.  (15))  takes  the 
form:  of  all  displacement  and  strain  functions  satisfying  prescribed 
displacements  on  Au  and  the  kinematic  Eqs.  (13),  the  displacement 
function  satisfying  the  equations  of  motion  (11),  the  constitutive 
equations  (12),  and  the  traction  boundary  conditions  on  AQ  are 
given  by 


6  {  f  rK  +  V2  +  »*<«?  +  . I  +  2^)  -  4  «£  *  « 

u  Jv 

"  [  L<nr  °r  +  ne  ^re)  ur  +  <°r  °r8  +  n9  50>  ueJ  "  0  <21) 


where 


n  =  nr  ^  +  n0  eQ 

U.  being  the  outward  normal  to  Aa  and  e^  ,  e^  the  unit  vectors  in 
the  radial  and  tangential  directions. 

Equation  (21)  is  obtained  from  Eq.  (15)  by  a  transformation  of 
coordinates  using  the  fact  that  and  are  invariants. 

The  alternative  form  of  Eq.  (21),  for  the  case  (for  simplicity) 
where  the  normal  n  -  (i.e.  for  a  circular  boundary)  is  obtained 

by  taking  variations  with  respect  to  ur  and  uq 


•fg 

“  (  [  ti?(X*(er  +  *e>  +  2*  er>  +  I  le^*  er0)  +  ^T(er  "  e6> 

Jv 

+  p  uj2  ur  j  6ur  +  &  (X*(er  +  e0)  +  1*  eQ)  +  2  |j(n*  «r0) 


* 

^  er.fl  9  *s  1 

+  — -  +  p  u)2  uj  6u0J 


dV  + 


(  [^<«r 


+  efi)  +  2p.  6 


-  a  "Uu 
r  rJ  i 


+  V*  er0  "  &r9J  5u9j  =  0 


(22) 


0  0 

where  and  ar9  are  the  prescribed  traction  components  on  . 

This  expression  is  seen  to  vanish  when  the  Eqs.  (11),  (12),  and 
the  traction  boundary  conditions  cn  Aa  are  satisfied. 

If  the  normal  points  inwards  (i.e.  ii  *  -e^)  the  sign  of  the  sur¬ 
face  integral  changes. 

Once  again,  if  the  trial  functions  for  displacement  components  are 
chosen  such  that  the  traction  boundary  conditions  are  satisfied  on  Aa 
the  surface  integral  vanishes  and  we  have  a  simpler  form  of  the  varia¬ 
tional  principle. 
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The  variational  principle  for  the  energy  equation  (18)  becomes 


{  f  [1  {(P)2  +7(I'/}  -  2  j  Ddl'ldV+j"  S  I  d*} 


where  D  is  given  explicitly  in  Eq.  (14). 

Its  alternative  form,  obtained  as  usual  by  taking  variations  with 
respect  to  the  temperature,  for  the  case  ®.  =  is 


-  f  (K  V*  T  +  2D)  6T  dV  +  f  (fi  +  K  |^)  «T  dA  *  0 

\ 


2 

where  V  is  as  given  in  Eq.  (14), 


Once  again,  as  expected,  this  expression  vanishes  when  the  energy 
equation  (14)  and  the  heat  flow  boundary  condition  on  Ajj  is  satis¬ 
fied. 


AN  EXAMPLE  IN  ONE  DIMENSION 

1.  The  Problem 

The  problem  of  steady  state  longitudinal  waves  in  a  viscoelastic 
rod  with  thermomechanical  coupling  is  now  solved  using  the  one  dimen¬ 
sional  versions  of  the  variational  principles  presented  in  the  previous 
section.  A  Rayleigh-Ritz  procedure  is  used  on  the  alternative  forms 
of  the  variational  principles.  The  same  problem,  including  time 
dependence,  was  solved  by  Huang  and  Lee  [9]  using  a  finite  difference 
approach.  The  results  obtained  here  are  compared  with  some  steady 
state  results  given  in  [9].  The  question  of  how  convergence  to  the 
solution  given  in  [9]  depends  on  the  number  of  coordinate  functions  used 


is  discussed. 


m  ^waw^Mggmwtyt  »com  jn»r  r» 


Let  us  consider  a  viscoelastic  rod  of  length  &  insulated  on  its 
lateral  surface  as  shown  in  Fig.  1.  The  left  end  is  free  while  the 
right  end  is  vibrated  at  a  frequency  ou  with  a  constant  stress  ampli¬ 
tude  cTq  (real),  so  that  the  prescribed  stress  at  this  end  is 
Oq  cos  <ut  .  The  temperature  of  the  vibrator  is  assumed  constant  at 
Tq  while  a  radiation  boundary  condition  is  assumed  at  x  =  0  . 

The  boundary  conditions  can  therefore  be  written  as 


x  *  0 


a  =  0 

§  -  «<*  •  V 


x  »  l 


a  *  a„ 


T  *  T„ 


where  c  »  h/K  is  the  ratio  of  the  surface  conductance  h  to  the 
thermal  conductivity  K  of  the  viscoelastic  material.  Note  that  the 
problems  of  uniform  normal  or  shear  traction  on  the  surface  of  a  wide 
slab  with  the  stated  thermal  boundary  conditions  prescribed  on  the  slab 
surface  are  mathematically  equivalent  problems.  Note  also  that  here 
we  have  mixed  thermal  boundary  conditions  but  this  can  be  taken  care 
of  in  the  variational  principle  as  shown  later. 


2.  The  Field  Equations  and  Variational  Principles 

The  equations  of  motion  (2)  and  the  constitutive  equations  (3) 
reduce  to 


^<E*^)  +  P«>2u  =  0 


CT  *  CTj  +  i  Og  »  E 


(27) 


-1 —  -  1  .  .1  1  - . 11  ------  - — — -  -»-  *.».  .•>.  »  ..  -■■■—  •■»«—■— ■»•  ■  ■»—'«««■ 
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where  E  *  E^  +  i  is  the  complex  Young’s  modulus  which  is  a  func¬ 
tion  of  the  temperature  through  the  reduced  frequency  (see  Eq.  (5)). 
The  steady  state  energy  equation  becomes 


K 


0 


(28) 


where,  as  before. 


du  du  du 

dx  dx  dx 


Note  that  for  the  one- dimensional  strain  problem,  E  must  be  re- 

•fc  Vfc 

placed  by  \  +  2(J-  and  E£  by  \ ^  +  2^  • 

The  corresponding  variational  principle  for  displacement  becomes: 
of  all  possible  displacement  functions,  the  one  satisfying  the  equa¬ 
tion  of  motion  (26),  the  constitutive  equation  (27)  and  the  stress 
boundary  conditions  from  Eq.  (25)  is  determined  from 


p  ui 


2-i 

u 


dx  +  cr 


u(*)j 


=  0 


(29) 


If  the  admissible  class  of  displacement  functions  is  restricted 
such  that  the  stress  boundary  conditions  are  already  satisfied  (through 
Eq.  (27)),  the  alternative  form  of  the  variational  principle  takes  the 
simplified  form 


M 

f  {  ^  VE*  S)  +  P  “2  UJ  6U  dX  "  C 


(30) 


The  temperature  variational  principle  takes  the  form:  of  all  tem¬ 
perature  distributions  which  satisfy  T(J&)  ■  Tq  ,  the  temperature  dis¬ 
tribution  which  also  satisfies  the  energy  equation  (28)  and  the 


section.  The  object  is  to  find  the  spatial  distribution  of  temperature, 
displacement  and  then  stress. 
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3.  The  Properties  of  the  Material 

A  Lockheed  solid  propellant  is  an  example  of  a  thermorheologically 
simple  materi;’.!  in  which,  within  a  wide  reduced  frequency  range,  the 
complex  shear  compliance  J  ■  -  i  J2  can  be  represented  by  (see 

[9]) 


"nl 

J1  =  ^ 

J2  «  k2(®  aT)  1 


where 


and  k^  ,  k2  ,  n^  ,  n2  ,  and  Tn2  are  constants.  The 
*  * 

compliance  D  ■»  1/E  is  related  to  the  shear  compliance 
* 

bulk  compliance  B  by  the  equation 


* 


D 


* 

3 


+ 


* 

B_ 

9 


tensile 

* 

J  and  the 


and  whenever  J  is  greater  than  B  by  at  least  two  orders  o Z  magni¬ 
tude,  we  can  write 


* 

D 


0  v 

°1  "  C1  "  Tl> 

B  Y 

°2  "  c2  (T  '  V 


(34) 


ttmtmrzzgz&T? 
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where  c^  ,  c^  ,  P  ,  Y  are  constants. 

E,  and  E,  are  now  obtained  from  above  as 

i  i 


E 


12  2 

DI  +  D2 


(35) 


'2  2  2 
Dl  +  D2 


4.  Method  of  Solution 

The  Rayleigh-Ritz  procedure  [11]  is  now  used  to  obtain  approximate 
solutions  for  the  temperature  and  displacement  (and  then  stress)  func¬ 
tions  from  the  variational  equations  (30)  and  (33). 

The  following  dimensionless  quantities  are  used 


T  -  T, 


q  --J 


T  = 


T  -  T  * 
A0  1 


(36) 


Equations  (30)  amd  (33)  are  nonlinear  and  it  is  not  possible  to 
choose  an  orthogonal  set  of  coordinate  functions  for  the  displacement 
and  temperature.  For  convenience,  it  is  assumed  that  the  displacement 
and  temperature  distributions  can  be  approximated  by  a  linear  combina¬ 
tion  of  polynomials  with  coefficients  to  be  determined.  Thece  func¬ 
tions  are  chosen  such  that  they  satisfy  the  boundary  conditions  (Eq. 
(25))  for  all  choices  of  these  unknown  coefficients.  Also,  the  number 
of  terms  in  the  series  are  parameters  which  can  be  set  in  the  result¬ 
ing  algebraic  equations  for  the  coefficients.  This  enables  comparison 
of  successive  approximations  with  the  solution  in  [9]  and  thus  an  esti¬ 
mate  of  the  degree  of  accuracy  as  a  function  of  the  number  of  terms 
taken  is  obtained. 


-  -  - - - 


The  non-dimensional  '  mperature  is  written  as 


a 


T(q)  -  1  +  (1  -  q)  {bQ  +  et  bQ  q  +  q2  £  b^l  -  q)1  2) 


(37) 


i=2 


where  eL  *  h  +  1  ,  m  is  a  parameter  and  bQ  »  b2  »  b3  "  *  bm  are 
m  real  constants  that  are  to  be  determined  (  m  <  2  implies 


Ft 

)  s  0  ). 

T*2 


It  is  easily  seen  thqt 


dT 

*>1,-0 


b0(e1  -  1)  -  H  b0  .  K(T  .  l)q=c 


and 


1  =  1 

U=i 

which  means  that  the  thermal  boundary  conditions  from  Eq.  (25)  are 
satisfied  in  terms  of  the  non-dimensional  variables  defined  in  Eq.  (36). 
The  complex  strain  is  written  as 


c  = 


C!  +  1  *2 


•%-*L 

i=0 


(38) 


where 


uQ  *  ag  +  i  aj  *  D*|  ^  ^  =  (ct  -  i  c2)  (TQ  -  T^  cQ 

n  is  a  parameter  and  ai  >  a2  *  a3  ’ *  *  an  are  n  comPlex  constants 
that  are  to  be  determined. 
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M  before,  It  Is  obvious  that  with  this  choice  of  strain 

o|  •  0  ,  "I  ■  E*|  ,  *0  “  °0 

!q«0  1 1*1  >  I'1 

so  that  the  stress  bonnier,  conditions  Iron  ft.  <»>  »  ““ 

Writing  a,  -  4  +  1  ^  C*  -  l,n)  this  choice  of  functions  leads 

to  (in  ♦  n>  real  unknowns  which  nust  be  determined  iron  an  e,»al  nu» 
ber  of  algebraic  equations. 

These  nonlinear  algebraic  equations  are  now  determined  free,  the 

. .  / qn\  and  (3 3\,  Substituting  the  displacement 

variational  equations  (30)  and 

,  p-  mm  cartvitiR  out.  the  necessary 
and  temperature  expressions  into  Eq.  (  )>  ^ 

integrations  and  equating  the  coefficients  of  »  a*  U  -  l.»> 
zero  gives 


for  j  “  1»2,3, . . .n 


2  —  +  )  a*  d[f(j,k)  -  f(j»k+D] 

(j+l)(j+2)(j+3)  ^  k 

k«0 

ft  a0  +  ^  S0>  ,T  .  ,T  .  +  I.  .  as.)  -  0  (39) 

-  i  —'  .  -72 -  (Ik+j  k+3+l  k+j+2 

U«sO  »  01 


-  a  -I  / a  —  i  n)  to  aero  gives 
and  equating  the  coefficients  of  aj  (j 


for  j  =  i»2,3, . . -n 


k«0 

cVj  -  + w 

k«0  * 30? 


•ownesetd^SHMCNi 
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where  Ip  (  p  an  integer)  is  a  nonlinear  function  of  e^  , 

b,  ,  b,  ,  b  defined  as 
z  j  tn 

1 

Ip  "  f  _ Q-l  q?P  dq _ _ _ 

0  i1  +  -  4bo +  ai  b0  * +  «2  L  V1  -  D1’2^ 

i=2 

f(j,k)  is  a  function  of  integers 


f(j,k)  = 


(  j+2)  (  J+3)  (k+1)  (k+2)  (k+j+3)  (k+j+4) 


d  is  a  non-dimensional  parameter 


2  2 
d  =  pj«ii! 


i  2  .  R.2  ,  .1.2 

|ao  s  <ao>  ■■  (a0> 


Note  that  Eq.  (30)  requires  the  displacement  u  in  addition  to 
the  strain.  Integration  of  Eq.  (38)  leads  to  an  extra  constant,  say 
Cq  ,  but  also  an  extra  equation  obtained  by  equating  the  coefficient 
of  6  Cq  to  zero.  This  extra  constant  Cq  has  been  eliminated  from 
Eqs.  (39)  and  (40)  given  above. 

Next,  substituting  the  displacement  and  temperature  expressions  into 
Eq.  (33)  and  equating  the  coefficient  of  6  bQ  to  zero  gives 


(3  +  et)  et  bQ 


»«  11 

1  bk^8k-l  "  2gk  +  8k+l*  +  V  I  L(ak  ai  +  ®k  ai* 

k*2  i,k**0 


{(l+e1) Il+k+i  ”(2  +  3at)  Ii+k+2+  +  3el>  Ii+k+3“  el  Ii+k+4^J  "  0 


and  equating  the  coefficients  of  6  b 


gives 


j  (j  *  2,m)  to  zero 


f°r  j-  2,3,4.., 


S  b0 


J(J+l)(J+2)  '  1,  Vh<J.k-l)  ‘  2tl(j ,k)  +  h(j,k+l)] 
k=2 

n 

+  V.^k*W  $  (Wl  -  *W  ♦  M 


i,k*=0 


i+j+k  i+j+k+1 


Ci+j+k+2  +  Ii+j+k+3i  *  0 


where  gk  and  h(j,k)  are  given  by 

(k+l4e  )(k-l) 

cr  *  — *• 

Bk  (k+1) 

h(J-k>  -7gr)-w-l)U-2) 

V  is  the  non-dimensional  parameter 

■■  ~*2  c2<T0  ~  »1'S 
2K(cJ  +  c*) 

and  Ip  has  been  defined  before  in  Eq.  (40). 

Equations  (39),  (40),  (41),  and  (42)  constitute  a  set  of  (2n-fm) 
nonlinear  algebraic  equations  for  the  (2n-hn)  unknowns  a£  (k  «  l,n)  , 
®k  (k  -  l,n)  ,  bQ  and  bR  (k  «  2,m)  . 

The  integral  Ip  given  in  Eq.  (40)  can  be  approximately  evaluated 
by  expanding  the  denominator  In  a  binomial  aeries  retaining  only  linear 


terms  and  carrying  out  the  integration. 


'OSHm sssa-v*^., 


For  example,  for  m  =  5  ,  we  get 


I  =  Y(1  +  V'1'"1  I 


wk  K!  p! 
(p+k+1)  ! 


where 


■<«!  "  b0 

-(b2  +  b3  +  b4  +  b5  -  eL  bQ) 
(b2  +  2b3  +  3b 4  +  4b5) 

-(b3  +  3b4  +  6b5) 

<\  + 


The  values  of  Ip  for  other  values  of  m  can  be  easily  calculated. 

This  approximation  for  Ip  proved  sufficiently  accurate  for  the 
calculations  carried  out.  I  can,  of  course,  be  more  accurately  deter¬ 
mined  by  numerical  integration  for  each  trial  value  of  b^  ,  b2  ,  . . . 
bm  during  iterative  solving  of  the  nonlinear  algebraic  equations  (39) , 
(40),  (41),  and  (42). 

The  stresses  are  determined  from  the  strains  and  temperature  from 


®1  *  E1  C1  -  E2  6 2 


°2  "  E2  S1  +  E1  €2 


and  the  stress  at  any  time 


-w  iwt 

Re(a(x,t))  *  Re(a  e  )  “  cos  “  °2  s*n 


These  equations  follow  inroediately  from  Eq.  (27).  and  E2  are 


determined  as  functions  of  temperature  from  Eqs.  (34)  and  (35), 
Non* dimensional  stresses  and  s2  are  defined  as 


and  at  x  ■  £ 


where 


81  “  ^  ai  *  s2  “  CT2 


s0  "  X  C0 


X  «  [2K  U)  p(TQ  -  Tl) ]" 


5.  Results  and  Conclusions 

Numerical  calculations  have  been  carried  out  for  the  following  data 
for  a  Lockheed  solid  propellant  [9]  whose  mechanical  and  thermal  proper¬ 
ties  are  qualitatively  typical  of  many  viscoelastic  solids 

cL  «  4.61  x  10'U  (psi)"1  (sec)P  (*F)"Y 
c2  *  1.62  x  10“U  (psi)"1  (sec)B  (°F)"Y 


P  -  -  0.214  y  -  3.21 

K  -  1.0  ^  Tq  *  65°F 

Tx  -  -  125°F  £  =  3  in. 

£2p  *  1.023  x  10  4  psi-sec2 

2Kp (Tq  -  Tj)  -  8.08  x  10"4  psi2- sec 

4 

w  »  10  rad/sec.  Sq  ■  0.5  (a^  «  1,42  psi) 


'in  [9]  h  should  read  1.0  instead  of  0.1. 
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The  nonlinear  algebraic  equations  (39) ,  (40) ,  (41) ,  and  (42)  were 
solved  numerically  in  a  computer  for  different  values  of  n  and  m  . 

The  subroutine  used  is  given  in  [14].  The  method  is  a  compromise  be¬ 
tween  the  Newton-Raphson  algorithm  and  the  method  of  steepest  descent. 

Figures  2,  3  and  4  show  the  resulting  t  ,  and  S2  distribu¬ 

tions  for  different  values  of  n  and  m  and  also  the  solution  from 
[9]  obtained  by  the  method  of  finite  differences.  The  solution  for 
n  *  1  ,  m  “  0  is  crude  but  we  see  that  the  convergence  to  the  true 
solution  is  very  rapid.  Figures  5,  6,  and  7  show  the  approximate  solu¬ 
tions  for  n  =  4  ,  m  =  3  .  Even  with  these  relatively  small  number  of 
terms,  the  stress  solutions  are  practically  identical  to  those  given  in 
[9],  while  the  temperature  solution  is  well  within  engineering  accuracy. 
The  algorithm  for  solving  the  nonlinear  algebraic  equations  converges 
very  quickly  and  more  accurate  solutions  can  be  obtained,  if  desired, 
by  taking  larger  values  of  n  and  m  . 

As  mentioned  earlier,  if  the  steady  state  values  of  stress  and 
temperature  are  of  interest  (this  is  often  the  case  in  design),  the  method 
used  here,  which  yields  the  steady  state  directly,  is  superior  to  that 
used  by  Huang  and  Lee  in  [9]  where  the  complete  time  histories  of  the 
above  mentioned  quantities  were  determined.  In  some  cases  in  [9]  the 
authors  obtained  the  steady  state  solutions  by  numerically  integrating 
forward  in  time  till  the  variables  of  interest  did  not  change  significant¬ 
ly.  In  other  cases,  they  did  not  integrate  upto  the  steady  state  but 
stopped  at  some  large  value  of  time. 


The  results  obtained  are  thus  very  satisfactory  as  long  as  s^  and 
S2  are  sufficiently  smooth  functions  so  that  approximation  by  a  series 
of  polynomials  is  efficient.  The  nature  of  the  spatial  distribution  of 
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stress  depends  upon  the  particular  choice  of  frequency  and  driving 
stress.  For  a  given  driving  frequency,  larger  driving  stresses  lead 
to  larger  temperatures  since  more  mechanical  energy  is  dissipated  as 
heat.  This  causes  the  material  to  become  softer,  so  that  lower  stress 
wave  velocities  and  therefore  lower  wave  lengths  result.  If  s^  and 
s^  are  rapidly  oscillating  functions  of  q  ,  the  polynomial  series  is 
no  longer  efficient  since  a  larger  number  of  terms  must  be  taken  to  get 
the  required  accuracy  and  the  lack  of  orthogonality  of  the  polynomials 
gives  rise  to  Hilbert  matrices.  This  results  in  convergence  problems 
for  the  algorithm  used  to  solve  the  algebraic  equations.  The  varia¬ 
tional  principles,  however,  should  work  fine  for  these  cases,  if,  for 
example,  trigonometric  functions  are  chosen  instead. 

To  sum  up,  the  variational  approach  seems  comparable  to  the  finite 
difference  approach  for  waves  in  one  dimension  and  ought  to  be  more 
efficient  in  two  or  three  dimensions  where  the  differential  equa¬ 
tions  are  partial  and  finite  difference  simulation  becomes  much  more 
complicated.  Solving  for  displacements  instead  of  stresses  has  the 
advantage  of  automatic  satisfaction  of  compatibility  conditions  and 
Mitchell's  equations  for  multiply  connected  regions.  The  choice  of 
coordinate  functions  is  very  important  and  must  be  made  carefully. 
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Figure  4.  S2(q)  curves  for  different  values  of 


Figure  5.  Comparison  of  T(q)  (n  =  4,  m  =  3)  with  correct  solution. 


